
clear
set more off
program drop _all

local owndir = "" // put your path here

global dirdata = "`owndir'/Raw_Data"
global dirin = "`owndir'/Dta"
global dirout = "`owndir'/Output"


************
* WEEKDAY EPSILON KAPPA NU ESTIMATION
************

*** 1. Run Regressions: Weekdays

use "$dirin/201911_W_clean.dta",clear
rename origin district
merge m:1 district using "$dirdata/district_list.dta"
drop _merge
drop district
rename d district_num
rename destination district
merge m:1 district using "$dirdata/district_list.dta"
drop _merge
drop district
rename d j
tsset district_num j
tsfill, full
bys district_num: egen total_number=sum(number)
gen pi=number/total_number*(1-0.17)
replace pi=0.17 if district_num==j // home sector (how can I add i-i commuting share?)
save "$dirin/temp.dta",replace

use "$dirdata/d_GIS.dta",clear
reshape long var, i(district_num) j(j)
merge 1:1 district_num j using "$dirin/temp.dta"
gen ln_pi=ln(pi)
gen d=var
reg ln_pi d i.district_num i.j
// -0.1413

*** 2. epsilon estimation with wage variation

use "$dirin/Seoul_Survey_clean.dta",clear

gen age_group=0
replace age_group=1 if de2>=2&de2<=5 // 20-30, 30-40, 40-50, 50-60
replace age_group=2 if de2==6 // 60+
drop if age_group==0 // 10-20

gen income=.
replace income=0.5 if de4==1 // < 1 million Won
replace income=1.5 if de4==2 // 1-2 million Won
replace income=2.5 if de4==3 // 2-3 million Won
replace income=3.5 if de4==4 // 3-4 million Won
replace income=4.5 if de4==5 // 4-5 million Won
replace income=6 if de4==6 // > 5 million Won

gen log_wage=ln(income)

egen destination_num=group(destination)
reg log_wage i.destination_num i.age_group
predict residual,r
gen exp_residual=exp(residual)

sum exp_residual

// in matlab calculate (0.4420/1.1019)^2 first (std.dev./mean)^2
// fminsearch(@(epsilon)abs(gamma(1-2/epsilon)/gamma(1-1/epsilon)^2-1-0.1609),5)
// estimated epsilon = 4.1642

// kappa=nu_weekday/epsilon_weekday=0.1413/4.1642=0.0339

************
* WEEKEND EPSILON KAPPA NU ESTIMATION
************

*** 1. Run Regressions: Weekends

use "$dirin/201911_H_clean.dta",clear
rename origin district
merge m:1 district using "$dirdata/district_list.dta"
drop _merge
drop district
rename d district_num
rename destination district
merge m:1 district using "$dirdata/district_list.dta"
drop _merge
drop district 
rename d j
tsset district_num j
tsfill, full
bys district_num: egen total_number=sum(number)
gen pi=number/total_number*(1-0.46)
replace pi=0.46 if district_num==j // home sector (how can I add i-i commuting share?)
save "$dirin/temp.dta",replace

use "$dirdata/d_GIS.dta",clear
reshape long var, i(district_num) j(j)
merge 1:1 district_num j using "$dirin/temp.dta"
gen ln_pi=ln(pi)
gen d=var
reg ln_pi d i.district_num i.j
// -0.1666

* epsilon_weekend=nu_weekend/kappa=0.1666/0.0339=4.9144

****
** FIVE PARAMETERS ESTIMATED HERE
****
// epsilon_weekday=4.1642
// epsilon_weekend=4.9144
// kappa=0.0339
// nu_weekday=0.1413
// nu_weekend=0.1666
